Exhausted CD8 T cells are associated with ICB responsiveness in breast cancer
Background
The single-cell dataset generated by Bassez et al. (2021) Nature Medicine contains single-cell data from breast carcinoma biopsies from 29 patients taken immediately before (pre-treatment) and 9 +/- 2 days after (on-treatment) receiving anti-PD1 immunotherapy. Lacking an objective measure of clinical benefit following therapy, in this study the authors use T cell expansion following treatment as a proxy for response.
Goal
To compare the T cell subtype composition between responders and
non-responders using ProjecTILs and reference maps for CD8 T cells and
CD4 T
cells.
R Environment
library(renv)
renv::restore()
library(patchwork)
library(ggplot2)
library(reshape2)
library(Seurat)
library(ProjecTILs)scRNA-seq data preparation
ddir <- "input/Bassez_breast_data"
if (!file.exists(ddir)) {
dir.create(ddir)
dataUrl <- "https://drive.switch.ch/index.php/s/7jdqnvPZuJrriZB/download"
download.file(dataUrl, paste0(ddir, "/tmp.zip"))
unzip(paste0(ddir, "/tmp.zip"), exdir = ddir)
file.remove(paste0(ddir, "/tmp.zip"))
}
# Count matrices
f1 <- sprintf("%s/1864-counts_tcell_cohort1.rds", ddir)
cohort1 <- readRDS(f1)
dim(cohort1)[1] 25288 53382
# Metadata
meta1 <- read.csv(sprintf("%s/1870-BIOKEY_metaData_tcells_cohort1_web.csv", ddir))
rownames(meta1) <- meta1$Cell
head(meta1) Cell nCount_RNA
BIOKEY_13_Pre_AAACCTGCAAGAAGAG-1 BIOKEY_13_Pre_AAACCTGCAAGAAGAG-1 1252
BIOKEY_13_Pre_AAACGGGTCGTAGGAG-1 BIOKEY_13_Pre_AAACGGGTCGTAGGAG-1 1869
BIOKEY_13_Pre_AAAGATGAGCAGGCTA-1 BIOKEY_13_Pre_AAAGATGAGCAGGCTA-1 1000
BIOKEY_13_Pre_AAAGATGCAAGCCGTC-1 BIOKEY_13_Pre_AAAGATGCAAGCCGTC-1 1288
BIOKEY_13_Pre_AAAGATGGTTTGACAC-1 BIOKEY_13_Pre_AAAGATGGTTTGACAC-1 2056
BIOKEY_13_Pre_AAAGATGTCAACGCTA-1 BIOKEY_13_Pre_AAAGATGTCAACGCTA-1 1224
nFeature_RNA patient_id timepoint expansion
BIOKEY_13_Pre_AAACCTGCAAGAAGAG-1 700 BIOKEY_13 Pre n/a
BIOKEY_13_Pre_AAACGGGTCGTAGGAG-1 995 BIOKEY_13 Pre n/a
BIOKEY_13_Pre_AAAGATGAGCAGGCTA-1 627 BIOKEY_13 Pre n/a
BIOKEY_13_Pre_AAAGATGCAAGCCGTC-1 681 BIOKEY_13 Pre n/a
BIOKEY_13_Pre_AAAGATGGTTTGACAC-1 789 BIOKEY_13 Pre n/a
BIOKEY_13_Pre_AAAGATGTCAACGCTA-1 634 BIOKEY_13 Pre n/a
BC_type cellSubType cohort
BIOKEY_13_Pre_AAACCTGCAAGAAGAG-1 HER2+ gdT treatment_naive
BIOKEY_13_Pre_AAACGGGTCGTAGGAG-1 HER2+ NK_REST treatment_naive
BIOKEY_13_Pre_AAAGATGAGCAGGCTA-1 HER2+ CD4_N treatment_naive
BIOKEY_13_Pre_AAAGATGCAAGCCGTC-1 HER2+ CD8_EM treatment_naive
BIOKEY_13_Pre_AAAGATGGTTTGACAC-1 HER2+ CD8_N treatment_naive
BIOKEY_13_Pre_AAAGATGTCAACGCTA-1 HER2+ CD8_RM treatment_naive
data.seurat <- CreateSeuratObject(cohort1, project = "Cohort1_IT", meta.data = meta1)Subset pre-treatment biopsies
table(data.seurat$timepoint)
On Pre
27528 25854
data.seurat <- subset(data.seurat, subset = timepoint == "Pre")Subset T cells according to annotation by the authors
table(data.seurat$cellSubType)
CD4_EM CD4_EX CD4_EX_Proliferating
4810 1707 107
CD4_N CD4_REG CD4_REG_Proliferating
3007 2861 69
CD8_EM CD8_EMRA CD8_EX
5443 290 2286
CD8_EX_Proliferating CD8_N CD8_RM
340 354 1963
gdT NK_CYTO NK_REST
1000 230 938
Vg9Vd2_gdT
449
data.seurat <- subset(data.seurat, subset = cellSubType %in% c("NK_REST", "Vg9Vd2_gdT",
"gdT", "NK_CYTO"), invert = T)We will remove samples that are too small and downsample large ones, in order to have similar contribution from all patients/samples. Downsampling large samples also speeds up downstream computations.
ds <- 1000
min.cells <- 200
tab <- table(data.seurat$patient_id)
keep <- names(tab)[tab > min.cells]
data.seurat <- subset(data.seurat, patient_id %in% keep)
Idents(data.seurat) <- "patient_id"
data.seurat <- subset(data.seurat, cells = WhichCells(data.seurat, downsample = ds))
table(data.seurat$patient_id)
BIOKEY_1 BIOKEY_10 BIOKEY_11 BIOKEY_12 BIOKEY_13 BIOKEY_14 BIOKEY_15 BIOKEY_16
1000 1000 555 1000 1000 939 1000 1000
BIOKEY_19 BIOKEY_2 BIOKEY_21 BIOKEY_24 BIOKEY_27 BIOKEY_28 BIOKEY_3 BIOKEY_31
1000 960 203 347 601 594 266 349
BIOKEY_4 BIOKEY_5 BIOKEY_6
1000 1000 616
ProjecTILs analysis
Load reference maps
Here we will project the query data onto human reference TIL atlas for CD4 and CD8 T cells, with the goal to interpret T cell diversity in the context of reference T cell subtypes. First, load the reference maps:
options(timeout = max(900, getOption("timeout")))
download.file("https://figshare.com/ndownloader/files/38921366", destfile = "CD8T_human_ref_v1.rds")
ref.cd8 <- load.reference.map("CD8T_human_ref_v1.rds")[1] "Loading Custom Reference Atlas..."
[1] "Loaded Custom Reference map Human CD8 TILs"
download.file("https://figshare.com/ndownloader/files/39012395", destfile = "CD4T_human_ref_v1.rds")
ref.cd4 <- load.reference.map("CD4T_human_ref_v1.rds")[1] "Loading Custom Reference Atlas..."
[1] "Loaded Custom Reference map custom_reference"
a <- DimPlot(ref.cd8, cols = ref.cd8@misc$atlas.palette, label = T) + theme(aspect.ratio = 1) +
ggtitle("CD8 T reference") + NoLegend()
b <- DimPlot(ref.cd4, cols = ref.cd4@misc$atlas.palette, label = T) + theme(aspect.ratio = 1) +
ggtitle("CD4 T reference") + NoLegend()
a | bReference-based analysis
We will project CD4 and CD8 T cells separately into the two reference
maps, then combine the results for a complete annotation of T cells
subtypes. For optimal batch-effect correction, we recommend projecting
each patient/batch separately (set
split.by="patient_id").
For large samples, Run.ProjecTILs can take several
minutes. Set the number of parallel jobs with ncores according
to your computational resources.
ncores = 8
cd4.projected <- Run.ProjecTILs(data.seurat, ref.cd4, ncores = ncores, split.by = "patient_id")cd8.projected <- Run.ProjecTILs(data.seurat, ref.cd8, ncores = ncores, split.by = "patient_id")See projected cells for two selected patients, projected onto the reference maps
s1 <- subset(cd4.projected, subset = patient_id == "BIOKEY_13")
a <- plot.projection(ref.cd4, query = s1, linesize = 0.5, pointsize = 0.2) + ggtitle("Patient 13 - CD4 T")
s2 <- subset(cd8.projected, subset = patient_id == "BIOKEY_13")
b <- plot.projection(ref.cd8, query = s2, linesize = 0.5, pointsize = 0.2) + ggtitle("Patient 13 - CD8 T")
s3 <- subset(cd4.projected, subset = patient_id == "BIOKEY_10")
c <- plot.projection(ref.cd4, query = s3, linesize = 0.5, pointsize = 0.2) + ggtitle("Patient 10 - CD4 T")
s4 <- subset(cd8.projected, subset = patient_id == "BIOKEY_10")
d <- plot.projection(ref.cd8, query = s4, linesize = 0.5, pointsize = 0.2) + ggtitle("Patient 10 - CD8 T")
(a | b)/(c | d)We can verify the expression profile using a panel of marker genes
genes4radar <- c("CD4", "CD8A", "TCF7", "CCR7", "IL7R", "LMNA", "GZMA", "GZMK", "FGFBP2",
"XCL1", "CRTAM", "TOX", "PDCD1", "HAVCR2", "PRF1", "GLNY", "GZMB", "TRAV1-2",
"KLRB1", "FOXP3")
plot.states.radar(ref = ref.cd8, query = cd8.projected, genes4radar = genes4radar)plot.states.radar(ref = ref.cd4, query = cd4.projected, genes4radar = genes4radar)Combine CD4 and CD8 T cell annotations in one set of cell type labels
data.seurat$functional.cluster.cd4 <- NA
data.seurat@meta.data[colnames(cd4.projected), "functional.cluster.cd4"] <- cd4.projected$functional.cluster
data.seurat$functional.cluster.cd8 <- NA
data.seurat@meta.data[colnames(cd8.projected), "functional.cluster.cd8"] <- cd8.projected$functional.cluster
annos <- data.seurat@meta.data[, c("functional.cluster.cd4", "functional.cluster.cd8")]
consensus <- apply(annos, 1, function(x) {
if (is.na(x[["functional.cluster.cd4"]]) & !is.na(x[["functional.cluster.cd8"]])) {
x[["functional.cluster.cd8"]]
} else if (is.na(x[["functional.cluster.cd8"]]) & !is.na(x[["functional.cluster.cd4"]])) {
x[["functional.cluster.cd4"]]
} else {
NA
}
})
data.seurat$functional.cluster <- consensus
data.seurat$functional.cluster <- factor(data.seurat$functional.cluster, levels = unique(data.seurat$functional.cluster))Compare subtype composition of responders vs. non-responders (based on clonal expansion) on pre-therapy samples
a <- ref.cd8@misc$atlas.palette
b <- ref.cd4@misc$atlas.palette
cols <- c(a, b)
by.exp <- SplitObject(data.seurat, split.by = "expansion")
tb <- table(factor(by.exp$E$functional.cluster, levels = names(cols)))
tb <- tb * 100/sum(tb)
tb.m <- reshape2::melt(tb)
colnames(tb.m) <- c("Cell_state", "Perc_cells")
p1 <- ggplot(tb.m, aes(x = Cell_state, y = Perc_cells, fill = Cell_state)) + geom_bar(stat = "identity") +
theme_bw() + scale_fill_manual(values = cols) + theme(axis.text.x = element_blank(),
legend.position = "none") + ggtitle("Responders") + ylim(0, 40)
tb <- table(factor(by.exp$NE$functional.cluster, levels = names(cols)))
tb <- tb * 100/sum(tb)
tb.m <- reshape2::melt(tb)
colnames(tb.m) <- c("Cell_state", "Perc_cells")
p2 <- ggplot(tb.m, aes(x = Cell_state, y = Perc_cells, fill = Cell_state)) + geom_bar(stat = "identity") +
theme_bw() + scale_fill_manual(values = cols) + theme(axis.text.x = element_blank(),
legend.position = "right") + ggtitle("Non-Responders") + ylim(0, 40)
p1 | p2The analysis shows that ICB-responding tumors have a higher proportion of CD8 exhausted (CD8.TEX), CD8 precursor exhausted (CD8.TPEX) and CD4 cytotoxic exhausted (CD4.CTL_Exh) T cells; on the other hand, CD4 Naive-like and CD8 TEMRA are more abundant in non-responders.
A more compact way of displaying enrichment/depletion of specific subtypes is to calculate and visualize fold-changes of TIL composition in responders vs non-responders
which.types <- table(data.seurat$functional.cluster) > 100 # disregard small populations
states_all <- names(cols)
cols_use <- cols[names(which.types)][which.types]
norm.c <- table(by.exp$NE$functional.cluster)/sum(table(by.exp$NE$functional.cluster))
norm.q <- table(by.exp$E$functional.cluster)/sum(table(by.exp$E$functional.cluster))
foldchange <- norm.q[which.types]/norm.c[which.types]
foldchange <- sort(foldchange, decreasing = T)
tb.m <- reshape2::melt(foldchange)
colnames(tb.m) <- c("Cell_state", "Fold_change")
pll <- list()
ggplot(tb.m, aes(x = Cell_state, y = Fold_change, fill = Cell_state)) + geom_bar(stat = "identity") +
scale_fill_manual(values = cols_use) + theme_bw() + geom_hline(yintercept = 1) +
scale_y_continuous(trans = "log2") + theme(axis.text.x = element_blank(), legend.position = "left") +
ggtitle("Responders vs. Non-responders")Tumors in ICB-responders have a >10 fold-change increase in CD4 exhausted CTL and a ~2-3 fold-change increase in both Tpex and **Tex*.
Correspondence analysis of subtype composition
Do patients cluster based on their subtype composition? are subtype composition “classes” related to response to ICB? we can perform correspondence analysis (CA) to try to answer these questions, and group patients based on their composition of T cell subtypes.
library("FactoMineR")
min.cells <- 30
# Only use samples with sufficient amount of cells
tab <- table(data.seurat$patient_id)
samples.use <- names(tab)[tab > min.cells]
data.seurat <- subset(data.seurat, subset = patient_id %in% samples.use)
subtype_comp_table <- prop.table(table(data.seurat$functional.cluster, data.seurat$patient_id),
margin = 2)
res.CA <- CA(t(subtype_comp_table), graph = FALSE)
a <- plot(res.CA) + theme(aspect.ratio = 1)
# Color by response
coords <- as.data.frame(res.CA$row$coord)
coords.ct <- as.data.frame(res.CA$col$coord)
coords.all <- rbind(coords[, 1:5], coords.ct[, 1:5])
is.R <- unique(data.seurat$patient_id[data.seurat$expansion == "E"])
is.NR <- unique(data.seurat$patient_id[data.seurat$expansion == "NE"])
coords[["Responder"]] <- "n/a"
coords$Responder[rownames(coords) %in% is.R] <- "Responder"
coords$Responder[rownames(coords) %in% is.NR] <- "Non-responder"
coords[["Patient"]] <- rownames(coords)
b <- ggplot(coords, aes(x = `Dim 1`, y = `Dim 2`, color = Responder)) + geom_point(size = 4) +
theme_bw() + theme(aspect.ratio = 1) + scale_color_manual(values = c("#E9E9E9",
"#fb0505", "#08e579")) + xlim(min(coords.all[, "Dim 1"]), max(coords.all[, "Dim 1"])) +
ylim(min(coords.all[, "Dim 2"]), max(coords.all[, "Dim 2"]))
a | bDiscussion
Our automated ProjecTILs analysis of this patient cohort indicate that compared to non-responders, early responders to PD-1 blockade (in terms of T cell expansion) have higher relative amounts of exhausted T cells (both CD8+ and CD4+) infiltrating their tumors, in agreement with previous studies (Daud et al. 2016; Thommen et al. 2018; Kumagai et al. 2020; Bassez et al. 2020). This supports the notion that patients with a pre-existing tumor-specific CD8 T cell response in their tumors are more likely to respond to PD-1 blockade therapy.
Because tumoricidal CD8 Tex are derived from Tpex (Siddiqui et al. 2019; Miller et al. 2019), the presence of CD8 Tex seems to be a good proxy for the presence of smaller populations of quiescent Tpex, that have the ability to proliferate and differentiate into Tex cells following PD-1 blockade (Siddiqui et al. 2019; Miller et al. 2019). Indeed we observed a small population of Tpex-like cells, with increased levels R vs RN. Because there are very few Tpex their gene profile might be less robust than others, but compared to Tex they display higher levels of Tpex-associated genes, such as XCL1 and CD200, and lower levels of terminal-exhaustion markers such as HAVCR2 and ENTPD1.
In their study, Bassez et al. also observed increased levels in R vs NR of a population they referred to as exhausted CD4 T cells. Here we observed also a trend for increased levels of follicular-helper CD4 T cells, that display exhaustion features (e.g. TOX, PDCD1, ENTPD1 expression), although the most prominent enrichment in ICB-responders is for exhausted CD8 T cells.
Further reading
Dataset original publication - Bassez et al
ProjecTILs case studies - INDEX - Repository
The ProjecTILs method Andreatta et. al (2021) Nat. Comm. and code
References
Bassez, A., Vos, H., Van Dyck, L. et al. A single-cell map of intratumoral changes during anti-PD1 treatment of patients with breast cancer. Nat Med 27, 820–832 (2021). https://doi.org/10.1038/s41591-021-01323-8
Daud, A.I., Loo, K., Pauli, M.L., Sanchez-Rodriguez, R., Sandoval, P.M., Taravati, K., Tsai, K., Nosrati, A., Nardo, L., Alvarado, M.D., Algazi, A.P., Pampaloni, M.H., Lobach, I.V., Hwang, J., Pierce, R.H., Gratz, I.K., Krummel, M.F., Rosenblum, M.D., 2016. Tumor immune profiling predicts response to anti-PD-1 therapy in human melanoma. J. Clin. Invest. 126, 3447–3452. https://doi.org/10.1172/JCI87324
Gros, A., Robbins, P.F., Yao, X., Li, Y.F., Turcotte, S., Tran, E., Wunderlich, J.R., Mixon, A., Farid, S., Dudley, M.E., Hanada, K., Almeida, J.R., Darko, S., Douek, D.C., Yang, J.C., Rosenberg, S. a, 2014. PD-1 identifies the patient-specific in filtrating human tumors. J. Clin. Invest. 124, 2246–59. https://doi.org/10.1172/JCI73639.2246
Miller, B.C., Sen, D.R., Al Abosy, R., Bi, K., Virkud, Y.V., LaFleur, M.W., Yates, K.B., Lako, A., Felt, K., Naik, G.S., Manos, M., Gjini, E., Kuchroo, J.R., Ishizuka, J.J., Collier, J.L., Griffin, G.K., Maleri, S., Comstock, D.E., Weiss, S.A., Brown, F.D., Panda, A., Zimmer, M.D., Manguso, R.T., Hodi, F.S., Rodig, S.J., Sharpe, A.H., Haining, W.N., 2019. Subsets of exhausted CD8+ T cells differentially mediate tumor control and respond to checkpoint blockade. Nat. Immunol. 20, 326–336. https://doi.org/10.1038/s41590-019-0312-6
Siddiqui, I., Schaeuble, K., Chennupati, V., Fuertes Marraco, S.A., Calderon-Copete, S., Pais Ferreira, D., Carmona, S.J., Scarpellino, L., Gfeller, D., Pradervand, S., Luther, S.A., Speiser, D.E., Held, W., 2019. Intratumoral Tcf1+PD-1+CD8+ T Cells with Stem-like Properties Promote Tumor Control in Response to Vaccination and Checkpoint Blockade Immunotherapy. Immunity 50, 195–211.e10. https://doi.org/10.1016/J.IMMUNI.2018.12.021
Thommen, D.S., Koelzer, V.H., Herzig, P., Roller, A., Trefny, M., Dimeloe, S., Kiialainen, A., Hanhart, J., Schill, C., Hess, C., Prince, S.S., Wiese, M., Lardinois, D., Ho, P.C., Klein, C., Karanikas, V., Mertz, K.D., Schumacher, T.N., Zippelius, A., 2018. A transcriptionally and functionally distinct pd-1 + cd8 + t cell pool with predictive potential in non-small-cell lung cancer treated with pd-1 blockade. Nat. Med. 24, 994. https://doi.org/10.1038/s41591-018-0057-z
Kumagai, S., Togashi, Y., Kamada, T. et al. The PD-1 expression balance between effector and regulatory T cells predicts the clinical efficacy of PD-1 blockade therapies. Nat Immunol 21, 1346–1358 (2020). https://doi.org/10.1038/s41590-020-0769-3